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ABSTRACT 


A  simplified  method  for  determining  the  aerodynamic  heating  that  occurs 
on  vehicles  traveling  at  high  speeds  is  presented.  Relations  are  submitted 
which  permit  the  calculation  of  the  temperature  of  the  boundary  layer  and  the 
heat  transfer  coefficient  as  a  <  (action  of  Mach  number  and  altitude.  Both 
laminar  and  turbulent  flow  withiu  the  boundary  layer  are  considered  and  the 
method  has  been  programmed  for  use  with  ABL  automatic  digital  computation 
equipment.  Equations  for  predicting  structural  temperatures  resulting  from 
aerod3mamic  heating  are  also  presented  and  are  incorporated  into  the  over-all 
program  in  order  that  the  effects  of  surface  temperature  on  the  aerodynamic 
heat  input  can  be  taken  into  account.  S 


The  program  can  also  be  applied  to  problems  other  than  those  specifically 
having  to  do  with  aerodynamic  heating.  Heat  flow  in  two  directions  accounting 
for  all  modes  of  heat  transfer  and  the  variation  of  physical  properties  V7ith 
temperature  may  be  considere<l^-~.AlsnL,_3tarious  materials  in  series  along  with 
resistances  at  their  interfacieAmay  be  studied  concurrently^^^ 


The  major  limitations  and  ratee  of  applicability  are  explained,  and  typical 
examples  are  presented  to  clarify  the  procedures.  ^ 


IMTRODUCTION 


A  vehicle  traveling  at  supersonic  speeds  within  the  atmosphere  generates 
severe  thermal  environments  by  its  high-speed  motion.  The  boundary  layer  concept 
proposed  by  Prandtl  has  proved  to  be  a  very  valuable  tool  in  studying  the  heating 
effects  due  to  flight  under  these  conditions.  The  particles  of  fluid  in  a  thin 
layer  near  the  body  are  acted  upon  by  shearing,  or  viscous,  frictional  forces 
and,  within  this  layer,  the  relative  velocity  between  the  vehicle  and  ambient 
fluid  is  slowed  doxjn  to  zero  at  the  surface.  The  kinetic  energy  of  the  air 
particles  then  appears  as  heat  energy  which  causes  the  boundary  layer  temperature 
to  rise.  Some  of  the  heat  from  the  boundary  layer  flows  into  the  surface  of  the 
vehicle,  the  amov^nt  increasing  rapidly  with  increasing  flight  speed.  This 
convective  heat-transfer  process  is  commonly  referred  to  as  aerodynamic  heating. 

The  obvious  result  of  aerodynamic  heating  is,  among  other  things,  a 
weakening  of  the  structure  if  it  is  not  adequately  protected  by  insulation  or 
some  means  of  cooling.  Aerodynamic  heating  is  of  particular  interest  to  the 
solid-propellant  rocket  designer  because  prolonged  exposure  to  such  heating  is 
likely  to  cause  a  failure  at  the  propellant -to-chamber  bond  surface  in  addition 
to  an  over-all  weakening  of  the  basic  load  carrying  structure.  It  is,  therefore, 
necessary  to  establish  a  means  of  evaluation  of  aerodynamic  heating  so  that 
proper  consideration  can  be  given  to  the  structural  design  and  heat  protection 
»  system  requirements  for  rocket  motors. 

The  purpose  of  this  report  is  to  describe  a  method  by  which  heat  transfer 
coefficients  and  boundary  layer  temperatures  can  be  calculated  and  in  turn  used 
^  for  determining  structural  response  to  aerodynamic  heating.  Both  laminar  and 

turbulent  flow  within  the  boundary  layer  have  been  considered  and  the  method  has 
been  programmed  for  use  V7ith  ABL  automatic  digital  computation  equipment. 
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The  procedure  is  an  intentionally  simplified  one  designed  to  provide  ansx^ers 
to  our  current  problems  and  should  not  be  expected  to  yield  satisfactory  results 
for  all  conditions  of  high-speed  flight.  Major  limitations  and  the  range  of 
applicability  are  explained.  Also,  typical  examples  are  presented  to  clarify 
the  procedures  involved  in  solving  several  types  of  aerodynamic  heating  and  heat 
transfer  problems. 

In  general,  the  method  can  be  used  for  flight  speeds  of  less  than  Mach  10 
and  altitudes  below  250,000  feet.  Very  high  boundary  layer  temperatures  result 
at  high  Mach  numbers;  consequently  the  oxygen  molecules  in  the  air  begin  to 
dissociate.  At  extremely  high  temperatures  (>.  10,000®R)  ionization  becomes 
significant.  Although  methods  are  available  for  treatment  of  these  phenomena, 
they  have  not  been  included  in  this  presentation. 

Interaction  of  the  boundary  layer  and  shock  wave  near  the  leading  edge  of 
a  hypersonic  vehicle  has  an  effect  on  the  aerodynamic  heating.  Consideration  of 
this  effect  was  not  considered  necessary  for  our  purposes.  Also,  radiation  to  or 
from  the  boundary  layer  has  been  neglected.  It  should  be  pointed  out  that  most 
of  our  (ABL)  studies  are  concerned  with  relatively  lovr,  compared  to  atmospheric 
re-entry  for  example,  velocity  flight  at  lcn7  altitudes  and  the  limitations  of  the 
present  program  need  not  concern  us  greatly.  For  example,  in  the  POLARIS  A3 
trajectory  the  results  may  appear  questionable  in  the  final  seconds  of  operation, 
due  to  the  limitations  previously  cited;  hov/ever,  the  peak  aerodynamic  heating 
and  maximum  skin  temperatures,  which  are  usually  the  items  of  major  concern,  are 
predicted  earlier  in  the  trajectory.  In  other  words,  by  the  time  such  effects  as 
dissociation  and  ionization  must  be  taken  into  account,  the  structure  is  losing 
energy  at  a  rate  greater  than  it  is  gaining  it.  Furthermore,  in  these  late 
stages  of  heating  when  the  coefficient  is  small,  structural  temperatures  are 
essentially  unaffected  by  its  changing. 
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The  surface  temperature  and  film  coefficient  of  heat  transfer  between  the 
boundary  layer  and  an  aerodynamically  heated  structure  are  very  closely  inter¬ 
related;  therefore,  it  was  necessary  to  include  a  procedure  for  determining  the 
temperature  distribution  within  the  structure.  Owing  to  its  inherent  simplicity 
and  speed  of  calculation,  the  two-dimensional,  finite-difference  technique 
presented  by  Dusinberre^^)  was  chosen  for  this  purpose.  Both  the  heat-transfer 
and  aerodynamic -heating  phases  of  the  resulting  program  can  be  used  independently; 
however,  only  an  approximate  value  of  the  film  coefficient  can  be  obtained  when 
the  surface  temperature  is  not  computed  at  the  same  time. 
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AERODYNAMIC  HEAi'  IMPUT 


Equations  representing  the  flow  of  a  compressible,  heat  conducting  and 
viscous  fluid  (such  as  air)  can  be  derived  from  the  principles  of  conservation 
of  energy,  momentum  and  mass,  together  with  an  equation  of  state.  These  equations 
can  be  simplified  by  introducing  the  boundary  layer  concept,  in  which  changes 
along  the  length  of  the  layer  are  considered  small  compared  to- the  changes  across 
it.  For  this  simple  case  the  equations,  based  on  energy  momentum  and  mass 
balances,  are: 


S-£li  +  =  0 

ax  ay 

The  equation  of  state  is: 


(1) 

(2) 

(3) 


p  RT  =  Constant  (3a) 

Consider  Equation  (1)  which  describes  the  conservation  of  energy.  On  the 

left  hand  side,  the  terms  in  parentheses  represent  convective  heat  transfer;  on 

the  right  hand  sid^,  the  first  temn  denotes  conductive  heat  transfer  and  the 

second  term  is  considered  a  "forcing  function"  which  is  proportional  to  au. 

ay 

For  a  small  boundary  layer  thickness, 6  ,  and  large  velocities,  u,  the  velocity 

gradient  within  the  layer,  au,  can  be  considered  approximately  equal  to  u. 

ay  S 


Therefore, 


A 

As  the  velocity  u  increases  /u\  assumes  an  increasingly  important  role.  This 

Ui 

forcing  function  represents  the  rate  of  dissipation  of  mechanical  energy  into 


heat  and  is,  therefore,  the  source  of  aerodynamic  heating. 
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The  above  equations  can  he  put  into  dimensionless  form  by  grouping  the 
variables  into  parameters  involving  two  or  more  variables.  The  advantage  of 
this  procedure  is  that  for  geometrically  similar  cases  the  equations  have  the 
same  values  for  the  parameters,  and  the  solutions,  in  dimensionless  form,  are 
equal . 

Equations  (1),  (2),  and  (3)  are  not  in  useful  form  since  many  of  the  terras 
are  difficult  to  evaluate.  A  major  effort  in  gas  dynamics  investigations  in 
recent  years  has  been  the  evaluation  of  derived  parameters,  such  as  the  familiar 
convective  film  coefficient,  for  suitably  representing  the  actual  process. 

Although  there  appears  to  be  almost  as  many  aerodynamic  heating  analysis 
methods  as  there  are  practicing  aerodynamic ists,  the  following  methods  have  been 
selected  and  developed  for  purposes  of  simplicity  and  engineering  application. 
These  procedures,  and,  consequently,  the  associated  computer  program  can  be 
modified  as  the  need  for  more  general  results  arises.  In  general,  results 
obtained  by  the  methods  reported  herein  will  be  valid  to  speeds  of  Mach  10  and 
altitudes  of  less  than  250,000  feet.  This  range  may  be  extended  by  incorporating 
relationships  to  account  for  dissociation  and  ionization. 

The  methods  and  equations  given  below  are  based  upon  the  "Reference 
Temperature  Method"  of  Eckert The  procedures  are  based  on  flow  past  a  flat 
plate  at  zero  angle  of  attack;  however,  the  same  relations  can  also  be  used  to 
approximate  conditions  on  slender  bodies  in  supersonic  flov?  whon  the  region  near 
the  leading  edge  is  excluded.  It  is  assumed  that  both  the  heat  transfer  and  skin 
friction  can  be  calculated  from  incompressible  flow  relations  provided  all 
temperature  dependent  air  properties  are  evaluated  at  an  appropriate  reference 
temperature  which  lies  somewhere  between  the  wall  and  the  boundary  layer  or 
recovery  temperature. 

Procedures  for  computing  the  boundary  layer  or  recovery  temperature,  the 
reference  temperature  and  the  film  coefficients  corresponding  to  both  laminar 
and  turbulent  flow  are  presented  in  the  following  subsections.  6 
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For  air,  T  =1.4,  and  Equation  (17)  reduces  to 

Tbi  =  Ibo  (1  +  0.17  MgZ) 


(18) 
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for  laminar  flow,  and 

(1  +  0.18  M^2)  (19) 

for  turbulent  flow. 

Heat  Transfer  Coefficients 

The  mechanism  of  heat  transfer  from  the  boundary  layer  to  a  surface  can  be 
written: 

q  =  h  A  (T„  -  (20) 

Thus,  in  estimating  the  heat  transfer  rate  to  a  structure,  it  is  necessary  to 
loim'r  the  heat  transfer  coefficient  beforehand. 

Laminar  Flow;  Assuming  constant  property  values  for  the  fluid,  it  has  been 
determined that  for  laminar  flow  across  a  flat  plate  v/ith  zero  pressure 
gradient,  the  local  skin  friction  factor  is 

Cf  =  0.664  Re"^  (21) 

and  the  local  Kusselt's  number  is 

Nu  =  0.332  (Re)%  (Pr)^^^  (22) 

By  substituting  Equations  (21)  and  (22)  into  Equation  (11),  the  Stanton  number 
for  laminar  flow  is  given  by 

St  =  2i  (Pr)"^'^  (22a) 

2 

The  Reynolds  analogy(^)  between  skin  friction  and  heat  transfer  is  expressed  as 

St  =  ^  (23) 

This  relationship  for  laminar  flow  is  a  reasonable  first  approximation  for 
subsonic  flow  since  the  Prandtl  number  is  approximately  1  and  it  doesn't  make 
a  great  deal  of  difference  where  the  fluid  properties  are  evaluated.  At  higher 
velocities,  however,  with  increasing  boundary  layer  temperatures,  the  Reynolds 
analogy  will  be  in  error  by  an  amount  approximately  equal  to  (Pr)"^/3  and  the 
fluid  property  values  become  increasingly  important. 
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For  laminar  flow,  the  heat-transfer  coefficient,  h,  can  now  be  written  as 


h  =  0.332  I  (Re)%(Pr)l/3 


(24) 


or 


h  =  0.332  (BS.)  ^  (Pr)^/^  (25) 

X 

For  supersonic  flight,  we  cannot  assume  constant  fluid  properties;  there¬ 
fore  in  evaluating  the  parameters  entering  into  the  solution  of  Equation  (25), 
the  variable  fluid  properties  must  be  taken  at  some  specific  temperature.  If  the 
temperature  chosen  is  the  boundary  layer  or  recovery  temperature,  we  find  that 
the  friction  factor  decreases  with  increasing  Mach  number  and  increasing  ratio 
of  .  On  the  other  hand,  if  the  temperature  chosen  is  the  wall  temperature, 

the  opposite  change  in  skin  friction  factor  is  noted.  From  this  it  might 
logically  be  inferred  that  some  intermediate  temperature  can  be  found  at  which 
the  variation  of  friction  factor  with  Mach  number  and  T^^/Too  vanishes.  It  is 
on  this  reasoning  that  the  Reference  Temperature  Method  is  based. 

Rubesion  and  Johnson(5)  have  shown  that  the  reference  temperature  can  be 
accurately  expressed  by  the  following  equation! 


T*  =  jj.  +  0.032  MgZ  +0.58 


*w 


\ 

-1) 


(26) 


This  expression  for  reference  temperature  can  be  rewritten  as 


T*  =  T.,-,  (0.42  +  0.032  M^^)  +  0.58  T^,  (27) 

The  last  factor  in  Equation  (27)  contains  the  term  T^  which  is  V7hat  we  ultimately 
want  to  determine.  However,  the  wall  temperature  T^  must  be  found  from  a  heat 
transfer  analysis  of  the  structure.  The  logical  way  to  approach  this  problem, 
and  the  one  used  in  this  instance,  is  to  divide  a  given  missile  trajectory  into 
increments  of  time  and  to  calculate  heat  transfer  and  surface  temperature  at  the 
end  of  a  time  increment  based  on  a  value  of  and  thus  h,  evaluated  at  a  surface 
temperature  at  the  beginning  of  the  time  increment.  The  newly  computed  value  of 
Tfj  is  then  used  to  determine  a  new  reference  temperature  and  film  coefficient  to 


be  used  over  the  next  increment  of  time 
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For  speeds  below  M5,  the  difference  between  boundary  layer  temperature  and 
wall  temperature  is  not  too  great,  so  with  small  error  we  can  substitute 
Equation  (18)  into  Equation  (27). 

T*  =  (0.42  +  0.032  +  o.58  [Tab  (1  +  0.17  llg^ 

and 

T*  =  Tqc.  (1  +  0.131  M^2)  (28) 

Thus,  in  certain  situations  we  are  able  to  determine  T*,  with  a  reasonable  degree 
of  accuracy,  independent  of  wall  temperature.  Equation  (28)  has  been  programmed 
for  use  in  calculating  the  film  coefficients  for  instances  when  the  problem  is 
such  that  the  wall  temperature  cannot  be  determined;  i.e.,  when  the  structural 
physical  properties  have  not  been  set,  such  as  in  preliminary  phases  of  rocket 
design.  When  possible,  the  more  desirable  method  of  determining  the  reference 
temperature,  Equation  (27),  is  used.  The  viscosity,  density  and  conductivity  of 
the  fluid  are  found  at  this  temperature  and  substituted  into  Equation  (25)  for 
computing  the  heat  transfer  coefficient  for  laminar  flow. 

Turbulent  F1ot7!  For  laminar  flow  the  calculation  of  friction  and  heat 
transfer  by  theoretical  means  is  usually  recognized  as  quite  reliable.  Such  is 
not  the  case  for  turbulent  flow,  however.  At  best,  turbulent  flow  theories  are 
semi -empirical  and  have  limited  ranges  of  validity.  Exact  solutions  of  the 
boundary  layer  equations,  based  on  the  present  day  understanding  of  the  problem, 
are  virtually  impossible.  This  necessitates  heavy  reliance  on  experimental 
measurements  to  obtain  friction  and  heat  transfer  data.  The  various  theories 
which  have  been  published  deviate  considerably  from  each  other,  especially  at 
higher  Mach  numbers.  The  following  procedure  agrees  V7ith  experimental  results 
within  the  investigated  range. 

For  a  flat  plate  at  zero  angle  of  attack  and  zero  pressure  gradient,  the 
following  formula  gives  good  agreement  with  experimental  data  up  to  Re  =  10^. 

^  =  0.0296  Re'°*^ 


(29) 
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Substituting  Equation  (29)  into  Equation  (22a)  for  the  Stanton  number 


St  =  =  0.0296  (Re)“0.2  (Pr)"2/3  (30) 

Therefore, 

Nu  =  Y  =  0.0296  (Re)0*8  (Pr)l/3  (31) 

Solving  for  the  film  coefficient 

h  =  0.0296  ^y)  (1*1)^^^  (32) 

As  in  the  case  of  laminar  flow,  the  fluid  properties  must  be  evaluated  at 
a  reference  temperature.  The  appropriate  reference  temperature  is  derived  in  a 
manner  analogous  to  that  for  laminar  flow.  Substituting  Equation  (19)  into 
Equation  (27), 

T*  =  T.^,.  (0.42  +  0.032  Ma2)  +  0.58  (1  +  0.18  Ma2)J  (33) 

which  simplifies  to 

T*  =  T  (1  4-  0.137  (34) 

Equations  (25)  and  (32)  indicate  the  form  of  solution  for  the  heat-transfer 
coefficient.  As  noted  earlier,  the  solutions  are  intended  to  be  applicable  to  a 
flat  plate  only.  For  laminar  flow  past  a  cone,  the  heat  transfer  computed  from 
Equations  (25)  and  (20)  should  be  multiplied  by  73.  In  the  case  of  turbulent 
flov7  past  a  cone,  multiply  Re  by  2  and  proceed  as  for  a  flat  plate. 

The  work  involved  in  picking  property  values  and  evaluating  the  quantities 
involved  in  the  solution  of  the  equations  for  the  film  coefficient  can  become 
most  laborious.  This  makes  hand  computations  over  a  trajectory  lasting  for  any 
appreciable  time  highly  impractical.  However,  the  problem  can  be  handled  very 
quickly  by  the  use  of  a  high-speed  digital  computer.  All  the  above  equations 
that  are  necessary  for  the  solution,  along  with  data  and  necessary  equations  for 
property  evaluation,  have  been  programmed  on  the  IBM  7074  computer  now  in  use  at 
ABL.  Details  pertaining  to  air  property  evaluation  and  computer  programming  are 
discussed  in  later  sections  of  this  report.  12 


Approximate  Method  for  Estimating  Heat  Transtcr  Coefficients 

Reference  (7)  contains  a  compilation  of  T^orking  curves  based,  with  some 
modification,  upon  the  reference  temperature  method  of  Eckert. This  method 
can  be  useful  in  making  computations  by  hand  if  an  accurate  estimate  of  the  wall 
temperature  can  be  made.  As  was  noted  earlier,  it  is  possible  in  the  lov7-speed 
regime  to  approximate  the  wall  temperature  with  the  boundary  layer  temperature. 

Another  very  simple  method  was  learned  by  the  writer  while  with  the  Martin 
Company  and  is  Included  for  very  rough  estimations  only.  An  equation  giving  the 
heat  transfer  index  has  been  prepared  In  nomograph  form  for  the  convenience  of 
the  designer  and  has  been  included  as  Figure  1.  The  following  procedures  are 
provided  to  explain  the  application  of  the  nomograph; 

1.  Plot  altitude  and  Mach  number  as  a  function  of  time  for  the  trajectory 
under  consideration. 

2.  Using  Standard  Atmosphere  Charts,  plot  ambient  temperature  as  a  function 
of  time. 

3.  Calculate  reference  temperature  from  Equation  (28)  and  plot  as  a 
function  of  time. 

4.  For  a  particular  time,  compute  the  ratio  of  reference  temperature  to 

ambient  temperature.  With  this  ratio  plus  the  Mach  number  and  altitude  at  that 

0  2 

time,  enter  Figure  1  to  determine  hx  *  where  x  is  the  distance  from  the  leading 
edge  or  nose.  Divide  the  index  by  x®*^  and  plot  on  a  curve  giving  h  as  a  function 
of  time.  The  h  values  can  then  be  used  to  make  a  crude  estimate  of  heat  transfer 
to  the  structure. 

Atmospheric  Properties 

Since  all  aerodjmamic  heating  calculations  are  dependent  upon  atmospheric 
temperature.  Figure  2  was  included  to  show  the  variation  of  ambient  temperature 
with  altitude  assuming  either  the  NACA  standard  atmosphere or  the  ANA  standard 
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'vjmo 'or  Approximate  CalculaHo.t  of  Aerodynamic  Heat  Transfer  Coefficient 


ALTITUDE,  THOUSANDS  OF  FEET 


FIGURE  2 

Varlaflon  of  Ambient  Temperature  with  Altitude 


THOUSANDS  OF  FEET 


hot  atmosphere.  To  simplify  the  computer  program,  the  temperature-altitude 
curves  were  conveniently  represented  by  algebraic  equations. 

The  variation  of  air  density  with  altitude  was  predicted  from  equations 
which  had  been  in  use  previously  in  related  studies at  ABL. 

Other  air  properties  necessary  for  calculating  the  aerodynamic  heat  transfer 
coefficient  are  thermal  conductivity,  viscosity  and  specific  heat.  Also,  since 
the  quantitative  values  for  each  of  these  properties  changes  with  changing 
temperature,  and,  since  they  must  be  determined  at  a  variable  reference  tempera¬ 
ture,  each  must  be  expressed  as  a  function  of  temperature.  This  also  holds  for 
the  air  density;  consequently  the  density  of  air  at  a  given  altitude  and  corres¬ 
ponding  ambient  temperature  must  be  modified  when  solving  for  the  density  at  the 
reference  temperature.  All  of  the  equations  relating  properties  and  temperature 
which  were  utilized  in  the  program,  along  with  the  reference  from  which  they 
V7ere  taken,  are  given  in  the  Equation  Summary  and  will  not  be  repeated  here. 

The  velocity  input  to  the  aerodynamic  heating  equations  must  be  in  terms  of 
Mach  number.  Since  the  velocity  in  feet  per  second  must  be  divided  by  the  speed 
of  sound  in  determining  Mach  number,  the  sonic  velocity  must  be  knoTO  for  the 
altitude  in  question.  The  atmospheric  temperature  is  available  from  Figure  2 
and  the  sonic  velocity  can  be  calculated  from  the  well  known  expression 

a  =  49.02  (T^3  )%  (35) 

The  sonic  velocity  up  to  an  altitude  of  100,000  feet  can  be  read  from  Figure  3. 


Ik 


ALTITUDI 


NUMERICAL  METHOD  FOR  TRANSIEOT  HEAT-FLOt'J  ANALYSIS 


The  general  numerical  method  of  Dusinberre  was  selected  as  the  most  simple 
and  convenient  procedure  for  determining  structural  response  to  aerodynamic 
heating.  The  development  of  the  fundamental  equations  is  briefly  described 
below.  Only  the  detail  considered  necessary  for  a  general  understanding  of  the 
problem  has  been  included  in  this  report;  those  wishing  a  more  thorough  treatment 
of  the  subject  are  referred  to  Reference  (1), 


Consider  the  point  0  as  shown  above.  The  temperature  at  that  point  is 
taken  as  representative  of  a  certain  region  which  includes  the  point,  and  the 
heat  flow  at  the  boundary  is  taken  as  the  summation  of  flows  from  surrounding 
regions.  The  transfer  can  take  place  by  any  of  the  three  modes  of  heat  transfer; 
that  is,  convection,  conduction  or  radiation.  In  the  general  case,  heat  can  also 
be  generated  within  the  region  itself;  however,  for  this  study  it  is  not 
considered  necessary  to  take  this  effect  into  account.  The  net  heat  flow  into 
the  region  must  be  stored  therein,  thus  causing  the  temperature  of  node  0  to  rise. 
Conduction 

For  a  conductive  interchange  from  some  node  i,  representative  of  a  region 
adjacent  to  the  region  containing  0,  the  rate  of  heat  transfer  qj^  to  node  0  is 
given  by 


‘Ik  “  ^  (*-0  “  '-i) 


(36) 
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where  is  the  conductance  along  the  path  i  to  0  and  t^  and  are  the  tempera¬ 
tures  at  nodes  0  and  i,  respectively. 

In  general,  for  a  single  material 

Kk  =  kS  (37) 

where  S  is  a  geometric  factor  depending  on  the  particular  subdivision  of  the 
system.  In  this  report  only  squares,  rectangles,  and  annular  networks  are 
considered;  the  S  values  corresponding  to  each  of  these  subdivisions  are  given 
in  Table  I.  The  program  is  not  limited  to  a  subdivision  consisting  of  the  above 
geometric  figures,  however.  For  situations  in  which  a  triangular  network  appears 
necessary,  the  geometric  factors  can  be  calculated  from  expressions  given  in 
Reference  (1),  The  procedure  may  also  be  obtained  from  the  author  upon  request. 

For  a  single  material  the  conductance  is  obtained  by  multiplying  the 
geometric  factor  by  the  material  thermal  conductivity;  however,  for  materials  in 
series,  air  gaps  and  contact  resistances  the  rules  for  addition  of  conductances 
must  be  applied.  These  procedures  are  explained  in  almost  any  good  text  on  heat 
transfer. 

Convection 

For  the  convective  mode  of  heat  transfer  from  i  to  0,  the  flow  rate,  q^  is 

qc  =  Kc  (to  -  ti)  (38) 

In  this  case  h  is  the  convective  film  coefficient  and  A  is  the  cross-sectional 
area  normal  to  the  flow  path.  All  other  terms  are  as  given  for  conduction. 

Unlike  the  thermal  conductivity  k,  it  will  be  necessary  to  have  h  vary  with 
time  and  surface  temperature.  Furthermore,  it  will  often  be  necessary  to  calcu¬ 
late  h  for  each  time  increment  as  the  wall  surface  temperature  changes.  Procedures 
for  calculating  h  in  a  special  application  are  given  under  the  section  of  this 
report  entitled  "Aerodynamic  Heat  Input." 

Radiation 

For  the  case  of  radiation  between  bodies 

^r  =  ^  Ve 


(39) 


TABLE  I 

CalculaHon  of  Geometric  Coefficient  S 


Note;  The  area  A,  normal  to  the  heat  flow  path.  Is  bounded  by  the  perpendicular 
bisectors  of  adjacent  heat  flow  paths.  The  geometric  coefficient  S  Is 
obtained  by  dividing  A  by  the  length  of  the  corresponding  path. 


Fg  Is  usually  called  a  configuration  factor  and  is  a  function  of  angles  of 
radiation,  area  and  distances  between  radiating  surfaces.  Fg  is  an  eraissivity 
factor  to  account  for  the  departure  of  the  two  surfaces  from  complete  blaclcness 
and  is  a  function  of  the  emissivities  as  well  as  the  configuration  of  the 
surfaces.  O'  is  the  Stefan-Boltzmann  constant. 

Cf  =  1730  X  10-12  Btu/(hr.'^sq.ft.y'’]^ 

The  temperatures,  represented  by  Tg  and  Tj  must  be  expressed  on  the  absolute 
scale . 

It  is  convenient  to  represent  radiation  heat  transfer  in  a  manner  similar 
to  convection  and  conduction  in  which  case 


=  Kr  (To  > 


hence, 


(40) 


Ky  =  h^A 

Xi7ith 

hr  =  O'  (41) 

The  configuration  factor  and  emissivity  factor  in  the  above  equations  were 
combined  into  a  single  term  for  computer  programming. 

Since  this  study  is  primarily  oriented  toward  aerodynamic  heating  analysis j 
it  will  usually  be  assumed  that  radiation  is  from  the  aerodynamically  heated 
surface  to  a  surface  at  0°  temperature.  This  is  the  usual  assumption  for  radiation 
from  a  body  traveling  through  space.  In  this  case 

q^  =  O'  ATg^  (42) 

or 


qr  =  h^A  (To  -  Ti)  (43) 

and 

hr  =  Cfe  T^3 

For  radiation  to  space,  will  always  equal  0  in  the  equation  for  q^. 
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Energy  Storage  and  Heat  Balance 

The  net  heat  flow  into  region  0  must  be  stored  in  that  region  causing  the 
temperature  of  node  0  to  rise 


‘Inet  ^  ®  =  ^stored  -  ^d  P  ^^o '  "  ^o^ 
vrhere  the  superscript  is  used  to  designate  the  temperature  at  the  end  of  a  time 

interval  A  0.  Other  quantities  are  identified  in  the  List  of  Symbols.  A  heat 
balance  on  the  cross  hatched  zone  shown  in  Figure  4  gives 


Co  (^o'  -  to) 
^io  =  AS 


(45) 


where  Cq  is  the  heat  capacity  defined,  for  a  unit  depth,  by  the  following 
equation; 


Co  =  V^p  Cp  (46) 

There  will  be  a  region  for  which  the  heat  capacity  will  be  assigned  to  each 
particular  node  0.  The  area  of  the  region  is  designed  as  V^j.  The  method  for 
computing  the  area  of  the  region  associated  with  nodes  in  rectangular,  square 
and  annular ly  segmented  subdivisions,  are  summarized  in  Table  II.  The  equations 
for  a  triangular  network  are  given  in  Reference  (1).  A  rule  to  remember  in 
solving  for  is  that  the  area  Vq  is  bounded  by  the  perpendicular  bisectors  of 
the  flow  paths  to  node  0.  Once  the  area  associated  with  a  particular  node  has 
been  determined,  the  heat  capacity  is  computed  from  Equation  46, 

t 

By  expanding  Equation  (45)  we  obtain 

^l,0^*'i  "  ^0^  ^2,0^^2  “  ^o^  '■"""”"^,o^^n  "  ^o^ 


22- 
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Equation  (47)  is  strictly  true  for  convection  and  conduction  only;  therefore,  to 
include  radiation  all  temperatures  will  be  expressed  on  the  absolute  scale. 
Thus, 

™  I  Kl.o  Ae  .  ^2.0  A®  .  Kn.o  A9  E  Kq  Ag 

To  =  * - —  +  a  ■  — )  (48) 
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TABLE  n 

Calculation  of  Heat  Capacity  Area  Vo 

Notes  Tlie  area  Vo  Is  bounded  by  the  perpendicular  bisectors  of  the  flow  paths 
from  the  node  under  consideration. 

♦ 


2973 
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where  K  represents  coefficients  for  nil  temperature  dependent  heat  flov^s.  As 
was  previously  pointed  out,  the  procedure  for  calculating  K  will  depend  on  the 
mode  of  heat  flow  along  the  particular  path  In  question. 

Criterion  for  Stability 

Equation  (48)  Is  the  basic  working  equation  for  the  solution  of  transient 
heat  flow  problems  by  the  numerical  technique  considered  herein.  By  applying  this 
equation  to  every  node  in  a  system,  it  is  possible  to  compute  the  temperature 
history  of  the  system  given  the  initial  conditions  and  given  the  boundary 
conditions  as  a  function  of  time. 

Tlie  characteristic  of  this  and  similar  "explicit"  techniques  v/hich  is  most 
often  pointed  out  by  those  critical  of  the  method  is  the  relatively  short  time 
intervals  A 9  which  must  be  used  to  insure  stability  of  the  calculations.  By 
close  examination  of  Equation  (48)  it  can  be  seen  that  for 

SKo  Ae 


(1  - 


-)  <0 


(49) 


Tg  depends  on  in  a  negative  sense.  This  is  impossible  physically;  therefore, 
the  generally  accepted  criterion  to  insure  stability  is 


^^o  < 

Co 


(49a) 


This  criterion  must  be  satisfied  for  every  node  in  the  system.  In  formulating  a 
computer  program  based  on  the  above  equations,  arrangements  should  be  made  for 
the  machine  to  make  a  check  of  stability  at  each  node  and  for  all  time  intervals. 
This  is  desirable  since  the  values  for  K  at  each  node  will  be  varying  with  time 
in  the  general  case. 

Division  of  a  System  into  Nodes 

The  term,  system,  as  used  in  this  report,  is  any  region  of  space  and  its 
contents  which  is  regarded  as  isolated  for  purposes  of  observation  and  analysis. 
The  object  is  to  determine  the  temperature  distribution  in  a  physical  system  in 
vjhich  heat  transfer  is  taking  place  by  calculating  temperatures  at  certain  pre¬ 
selected,  discrete  points  or  nodes.  The  first  step  in  the  numerical  analysis  is 
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to  divide  the  system  by  a  network  of  reference  points  or  nodes  with  points  at 
which  temperatures  are  desired  included  eunong  the  reference  poitits.  Obviously, 
additional  nodes  may  also  be  necessary  particularly  in  areas  where  gradients  are 
changing  sharply.  As  was  described  above,  the  division  of  time  and  the  system 
subdivision  are  related  through  the  stability  criterion. 

The  method  described  in  this  report  can  be  applied  to  any  geometric  arrange¬ 
ment  of  nodes  as  long  as  the  region  assigned  to  each  node  and  the  heat  flov;  paths 
to  and  from  each  node  can  be  adequately  described.  Procedures  for  handling 
systems  subdivided  by  a  network  of  squares,  rectangles,  annular  segments  or  any 
combination  of  these  are  given  above;  triangular  networks  are  described  in 
Reference  (1). 

A  system  should  be  divided  as  simply  as  possible  using  as  few  nodes  as 
possible  since  each  additional  point  adds  to  the  time  and  labor  of  computation. 

The  computer  input  information  is  usually  simplified  if  a  given  subdivision  can  be 
repeated  throughout  a  large  portion  of  the  system  since  the  geometric  character¬ 
istics  need  only  be  entered  once  in  the  input  to  the  computer. 

It  was  mentioned  previously  that  the  usual  rules  for  addition  of  conductances 
can  be  applied  to  different  materials  in  series.  Although  this  is  true,  a  much 
simplier  method  can  be  utilized  if  a  judicious  selection  of  the  node  positions 
is  made.  Consider  two  materials,  1  and  2,  each  having  different  material 
properties  and  having  a  common  interface.  Also,  consider  the  heat  flow  from  a 
node  in  material  1  at  a  distance  from  the  interface,  to  a  node  in  material  2 
located  at  a  position  y2  from  the  interface.  Using  the  subscripts  1  and  2  to 
designate  materials  and  (f)  to  represent  the  conductance  across  the  interface, 
the  heat  flov7  by  conduction  from  node  1  to  node  2  is 

=  y^^ - 7~yJ  (T2  -  Tj)  (50) 

kj  f  +  ^ 
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and  the  conductance  is 


Ki,2  = 


yi  .  1  2  y2 
h  ^  h 


(51) 


This  value  must  be  calculated  and  entered  as  input  to  the  computer. 

The  simplier  method  is  to  position  the  nodes  at  the  interface  with  node  1 
considered  as  being  in  material  1  and  node  2  in  material  2.  Heat  transfer 
between  the  nodes  is  assumed  to  take  place  only  by  convection  across  the  inter¬ 
face  (this  is  reasonable  since  the  interface  can  be  assumed  to  absorb  no  heat) , 
In  this  arrangement  the  simplier  and  straight  fon^ard  methods  of  calculating 
conductances  can  be  applied. 
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SUMMARY  OF  EQUATIONS  AND  RELATED  IMFORtdATlON 


The  following  equations,  which  enter  into  the  combined  heat -transfer  and 
aerodynamic  heating  program  previously  described,  have  been  programmed  for  the 
IBM  7074  computer  at  ABL. 

Aerodynamic  Heat  Input  Reference  (2) 

1.  Boundary  Layer  Temperature 

(a)  Laminar  Flow 

a  +  0-17  Ma^) 

(b)  Turbulent  Flow 

Tfei  =  Toe,  (1  +  0-18 

2 .  Film  Coefficient 

(a)  Laminar  Flow 

h  =  0.332  (Re*)^(p/)^^^ 

(b)  Turbulent  Flow 

h  =  0.0296  ~  (Rg*)®'® 

Note:  The  character*,  appearing  as  a  superscript  over  a  given  property  indicates 
that  that  particular  property  is  to  be  evaluated  at  the  reference 
temperature, 

3.  Reference  Temperature 

(a)  General 

T*  =  T,;,^  (0.42  +  0.032  M^2)  +  o.58  T^^ 

(b)  Approximate  -  Laminar  Flow 

T*  =  (1  +  0.131  Mg2) 

(c)  Approximate  -  Turbulent  Flow 

T*  =  Too  +  0*T37 
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4,  Atmospheric  Properties 

(a)  Thermal  Conductivity 

k*  =  0.753  X  10"^  +  6.319  x  10"5  T* 

(b)  Viscosity 

_7  /  7 

n*  =  0.231  X  10-7  !  j 

(c)  Specific  Heat  at  Constant  Pressure 


Cp*  =  32.17  r 0.367  -  500  j 

T*  +  3500 

(d)  Velocity  of  Sound 

a  =  49.02  (T^ 

(e)  Density 


where 


for 


for 


=  p  ^  (1  -  6.89  X  10“^)^'^^® 


H  'k  35,352  ft.,  and 


_ 1.32 

exp  1 .432  + 


±o _ 

H-35.352 

20,950 


H>  35,352  ft. 


5.  Dimensionless  Numbers 


(a)  Reynolds  Number 


Re* 


P*  V 
^Ji  * 


(b)  Mach  Number 

Vs 

M  =  ~ 

^  ®s 

(c)  Prandtl  Number 

C„* 


k* 


Reference  (10) 


Reference  (11) 


Reference  (12) 


Reference  (13) 


Reference  (9) 


Reference  (9) 
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P  *  = 
r 


6 .  Atmospheric  Temperature 


Several  equations  X7ere  used  to  represent  the  temperature  as  a  function  of 
time  for  both  the  Standard  and  Standard  Hot  atmospheres.  These  are  given  in  the 
Abstract  of  the  Computer  program. 


Heat  Transfer 

1.  Node,  Temperature 

*^0'  =  +  ^2,0^2  + 


where 


and 


^i,o  ”  c 


AS 


0,0 


=  1  -  ^6 


Reference  (1) 


^n,o^n  ^0,0^0 


2 .  Conductance  Factors 
(a)  Conduction 
K  =  kS 


Note;  The  value  S  is  a  geometric  factor  depending  upon  the  geometric  division 
into  nodes  and  must  be  entered  as  input  to  the  computer.  The  area  A, 
normal  to  the  heat  flov?  path,  is  bounded  by  perpendicular  bisectors  of 
adjacent  paths.  The  factor  S  is  found  by  dividing  A  by  the  length  of  the 
corresponding  path. 

(b)  Convection  or  Radiation 
K  =  hA 


3.  Heat  Capacity 

Co  =  P  Cp 

Note;  The  area  Vg  is  bounded  by  the  perpendicular  bisectors  of  the  flow  paths 
leading  from  the  node  under  consideration.  This  term  must  be  submitted 
as  input  to  the  machine. 
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4.  Net  Heat  Transfer  to  the  Structure 

(a)  Rate  of  Transfer  (Aerodynamic  Heating) 
q(Qj  =  h(T|jj^  -  T„)  -  hj.(T„  -  T^) 

<b)  Total 

^  ^  0  ^(6) 

5.  Stability  Check 

Computer  automatically  adjusts  ^  0  until 

Z  ^  Q  g  fpj.  nodes  and  time  increments. 

Co 

Thermal  Properties 

Materials  thermal  properties  necessary  for  computing  transient  heat  flow  and 
structural  response  are  input  terms  and  can  be  entered  into  the  program  as  a 
function  of  temperature.  Some  typical  thermal  property  values  for  some  commonly 
used  Insulation  materials  are  given  in  Table  III,  These  are  considered  as 
suitable  for  use  in  preliminary  design  and  study  work;  however,  for  project 
design  work  and  more  detailed  analyses  involving  extended  temperature  ranges, 
these  values  are  not  recommended. 

Table  IV  is  a  listing  of  the  latest  thermal  property  data  for  Spiralloy 
which  are  known  at  the  present  time.  Since  work  is  being  continued  on  the 
evaluation  of  this  material,  and  since  the  properties  obviously  must  depend  on 
changing  fabrication  techniques  and  component  materials,  they  are  likely  to 
change  quite  rapidly.  In  the  past  the  Polaris  Division  has  been  the  leader  in 
obtaining  Spiralloy  thermal  property  data;  it  is  suggested  that  they  be  contacted 
when  a  detailed  thermal  analysis  of  a  Spiralloy  structure  is  required. 

Units 

Many  of  the  vjorking  equations  involve  constants  V7hich  were  determined  for  a 
given  set  of  units.  For  example,  to  account  for  radiation  all  temperatures  were 
expressed  in  degrees  Rankine,  Most  of  the  terms  entering  heat  transfer  calcu¬ 
lations  have  derived  units;  e.g,,  units  for  thermal  conductivity  could  be 
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_ Btu _  ,  and  can  be  accurately  expressed  using  several  combinations 

(sec,)(sq.ft,)®F/ft, 

of  the  fundamental  units.  Therefore,  to  avoid  confusion  and  to  avoid  the 
programming  of  numerous  conversion  factors,  the  equations  and  program  in  general 
have  been  set  up  for  the  so-called  English  system  of  units  using  Btu's,  feet,  and 
seconds  as  the  units  for  heat,  length  and  time,  respectively.  All  input  should 
be  reduced  to  this  one  consistent  system  of  units  for  calculation.  The  basic 
units,  along  with  the  derived  units,  are  identifiable  in  the  List  of  Symbols. 

When  transient  heat  transfer  calculations  are  carried  out  independent  of  the 
aerodynamic  heating  equations  and  radiation  effects,  the  above  procedure  concern¬ 
ing  units  need  not  necessarily  be  followed.  In  this  case,  all  that  is  required 
is  that  a  consistent  set  of  units  be  used. 
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Table  III 


Room  Temperature  Thermal  Properties  of  Some  Commonly 
_ Used  Insulation  Materials _ 


Thermal  Conductivity 
Btu 

Density 

lb. 

Specific  Heat 
Btu 

Material 

(sec.)(sq,ft.)(°P/ft.) 

■in  WHIM 

Reference 

Buna-S  Rubber 

7  X  10“5 

85.5 

0.44 

19 

Thermo -lag 

2.22  X  10“5 

SO.O 

0.298 

20 

Cellulose  Acetate 

3.8  X  10"5 

81.0 

0.32 

21 

Cork  (2755*) 

1.46  X  10‘5  ** 

35.6 

0,44 

18 

Boric-Acid  Filled 
NBR  Phenolic 

4  X  10-5 

77.4 

0.40 

18 

ATJ  Graphite 

1.458  X  10"2 

105,5 

0.243 

22  &  18 

RPD  -  150 

7.15  X  10-5 

109.2 

0.29 

23 

Pyrolytic  Graphite 

5.55  X  lO'J  (A-Plane) 
3.33  X  10*^  (B-Plane) 

125.0 

0.232 

24 

*  Armstrong  designation 

**  The  thermal  conductivity  of  this  type  of  cork  falls  quite  rapidly  up  to  a 
temperature  of  approximately  180°  at  which  k  =  1,076  x  10"^. 
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Table  IV 


Thermal  Properties  of  Splralloy 


Temperature 

°F 


Thermal  Conductivity^^®^ 
Btu 

(sec.)(sq.£t.)(°F) 


Density 


lbs, 

.  (cu.ft.) 


Specific  Heat (13) 
Btu 

(lbO(°F)- 


80 

7.70 

120 

7.92 

160 

8.14 

200 

8.37 

240 

8.61 

280 

8.88 

320 

9.12 

124 

0.248 

* 

* 

*  Data  not  available  as  a  function  of  temperature. 
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APPENDIX 


ILLUSTRATIVE  EXAMPLES 


The  following  examples  have  been  included  to  illustrate  the  method  of 
solution  of  heat  transfer  problems  by  the  procedures  presented  above  and  to 
point  out  the  capabilities  of  the  computer  program.  The  actual  preparation  of 
data  for  machine  input  is  explained  in  more  detail  in  Reference  14.  In  presenting 
the  results  for  each  example,  an  attempt  is  made  to  point  out  particular  program 
characteristics  relating  to  that  particular  type  of  problem. 

Example  I  is  a  typical  cooling  problem  for  which  a  ’’semi-analytical'' 
solution  is  readily  available.  The  results  from  the  computer  program  are 
compared  graphically  with  the  solution  taken  from  charts  which  correspond  to 
the  analytic  solution  and  which  are  presented  in  the  literature A  favorable 
comparison  is  seen  to  result  for  this  particular  problem.  Since  the  intent 
of  this  program  was  not  to  make  a  comprehensive  study  of  the  various  aspects  of 
heat  transfer  but  to  provide  a  simple  tool  for  aerodynamic  heating  analysis,  a 
more  detailed  study  and  further  comparisons  do  not  appear  necessary  at  this  time. 

Example  II  is  representative  of  several  aerodynamic  heating  problems  for 
vjhich  this  program  is  ideally  suited.  In  addition  to  pointing  out  the  type  of 
information  available  from  the  computer,  this  example  also  serves  to  illustrate 
the  effects  of  changing  cross-sections  and  the  importance  of  a  two-dimensional 
analysis  in  this  type  problem.  Ttro  approximate  methods  for  determining  film 
coefficient  are  also  compared  with  the  suggested  method. 

Example  III  was  included  to  show  the  comparison  of  computed  temperatures 
with  those  obtained  from  an  actual  flight  test.  Since  most  of  the  information 
pertaining  to  this  flight  is  classified,  only  the  temperature  results  are 
presented. 
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EXAMPLE  1 


Cooling  of  an  Infinite  Slab 

A  large  slab  of  unalloyed  titanium  4.25  inches  thick  is  held  at  1500®F  before 
being  hung  vertically  in  70®  air  to  cool.  A  fan  forces  air  over  the  slab  at  a 
velocity  such  that  the  convective  film  coefficient  is  maintained  at  a  constant 
62,7  Btu/(sq.ft.)(hr.)(®F),  It  is  desirable  to  know  the  surface  and  midplane 
temperature  as  a  function  of  time  until  the  surface  has  dropped  below  500®F. 

(a)  Consider  average  physical  properties  and  no  heat  loss  due  to 
radiation. 

(b)  Compute  new  temperature-time  relationships  if  radiation  is 
considered.  (Let  6  =  F^  x  Fg  =  0.7.) 

(c)  Determine  temperature  vs.  time  considering  variation  of  thermal 
properties  and  accounting  for  radiation. 

Solution 

Since  heat  is  being  transferred  at  equal  rates  from  both  sides  of  the  slab, 
only  the  temperatures  on  one  side  of  the  midplane  are  computed.  The  slab  is 
divided  into  9  elements  as  shown  in  the  following  sketch: 
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Each  element  Is  1/48  ft.  thick  by  1  ft.  square  except  for  the  element 
adjacent  to  the  surface  v/hich  is  1/96  ft.  thick.  This  results  in  a  total  of 
9  nodes  in  the  metal  plus  1  in  the  boundary  layer.  These  are  numbered  consecu¬ 
tively  from  the  outside  inward  starting  at  node  1  in  the  boundary  layer.  Node  10, 
although  not  exactly  at  the  midplane,  is  assumed  to  represent  midplane  temperature. 

A  listing  of  the  pertinent  property  values  is  given  below: 


Temperature 

Thermal 

Conductivitv 

Specific 

Heat 

80 

10.1 

0.125 

200 

9.8 

0.130 

400 

9.72 

0.136 

600 

9.82 

0.143 

800 

10.0 

0.151 

1000 

10.5 

0.160 

1200 

11.1 

0.171 

1400 

11.9 

0.183 

1510 

12.4 

0.19 

The  density  P  s  282  was  assumed  constant  throughout  the  range  of 


temperatures. 


The  thermal  properties  for  solving  parts  (a)  and  (b)  of  this  example  were 
taken  at  an  average  temperature  of  approximately  1200°F. 


The  follo\Jing  data  were  compiled  for  computer  input: 


1-  fin  n  i>UU 

(hr.) (sq.ft .T(^ 

(hr.)(sq.ft.)(“F/ft.) 
p  =  282  Ibs/cu.ft. 

•  '>•*58  (IbOCF) 

A  e  =  0.002  hr. 


Tbi  =  70°F 


Tj  (initial  temperature,  all  nodes)  = 
ISOO^F 

e  =  F^  X  F  =  0.7  (for  radiation 
part  of  the  problem) 


9  =  0.16  hr.  (This  value  was  taken  from  the  analytical  solution;  hov;ever,  the 
computer  can  be  made  to  stop  when  a  particular  node  has  reached 
a  certain  temperature,  if  desired.) 


Print  out  frequency  -  Arbitrarily  selected  as  every  tenth  time  increment, 

Vq  =  Nodal  volume 

V2  =  (l)(l)(l/96) 

Vq  for  nodes  3-10  is  0,02083  cu,ft. 

S  =  Cross-sectional  area  of  the  path  between  nodes  divided  by  the 
path  length, 

Sj^2  “  Mode  of  heat  transfer  is  convection;  cross-sectional  area 
(Aj^2  =  1)  required, 

S22  and  all  others  -  Mode  of  heat  transfer  is  conduction  and  S  =48  ft. 

The  heat  transfer  program  was  set  up  to  conserve  conqjuter  storage  space  by 
making  possible  the  division  of  a  system  into  five  separate  regions,  each  of 
tirhich  could  be  further  divided  by  the  nodal  network  into  a  maximum  of  50  nodes. 

As  a  result,  it  is  possible  to  enter  data  for  an  entire  region  without  having  to 
submit  repetitious  data  for  all  the  nodes  In  that  region.  In  general,  therefore, 
it  is  wise  to  lump  all  adjacent  nodes  having  similar  input  information  into  a 
single  region.  Each  region  is  ordinarily  connected  to  adjacent  regions  by 
conduction  between  nodes  at  their  common  boundary. 

Another  secondary  advantage  arose  out  of  the  above  method  of  programming. 

By  not  connecting  the  regions  in  a  given  run,  it  becomes  possible  to  run  five 
different  problems  at  one  time,  submitting  only  the  data  for  each  problem  (or 
region)  which  is  different  from  the  other  regions.  In  the  example  nov;  under 
consideration,  parts  (a),  (b),  and  (c)  were  all  run  at  one  time  and  took  approxi¬ 
mately  1  minute  including  print -out. 

The  results  for  Example  I  are  given  in  the  figures  which  follow.  Figure 
A-1  shows  the  surface  temperatures  as  a  function  of  time  for  the  three  different 
sets  of  conditions  which  were  specified.  Figure  A-2  presents  the  same  information 
but  for  the  midplane  temperature.  For  each  of  the  two  positions  and  for  average 
material  properties,  temperatures  were  also  computed  from  Shack  Charts which 
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are  based  on  analytical  solution  of  the  governing  differential  equations.  The 
results  are  shown  plotted  on  Figures  A-1  and  A-2,  and  are  seen  to  be  in  quite 
close  agreement  with  the  numerical  solutions. 

In  preparing  input  for  this  problem,  a  time  increment  of  0.002  hours,  known 
to  be  unstable  for  the  particular  subdivision  used,  V7as  selected.  The  computer 
automatically  reduced  the  time  increment  by  cutting  it  in  half  until  a  stable 
value  (  Aq  =  0.0008  hours)  was  reached. 


EXAMPLE  II 


Aerodynamic  Heating  of  a  Hypothetical  Air-to-Air  Missile 

Figure  A-3  shot/s  the  Mach  number  and  altitude  versus  flight  time  for  a 
hypothetical  air-to-air  missile.  The  trajectory  shown  is  assumed  to  follow  15 
minutes  of  flight  at  Mach  equals  1.47  at  40,000  feet  on  a  hot -atmosphere  day. 

The  following  information  is  to  be  computed  for  the  propulsion  unit  which  is 
shoTm  schematically  in  Figure  A-4. 

(1)  h  vs.  0  at  X  =  5.5  using  the  approximate  method  is  assumed 
equal  to  T^p. 

(2)  h  vs.  0  at  X  =  3.5,  5.5  and  7.5  taking  into  account  the  computed 
wall  temperatures. 

(3)  Tj^j^  vs.  0 

(4)  Propellant -insulator  interface  temperature  vs,  time  at  x  =  3,5, 

5.5  and  7,5. 

(5)  q^gj.  vs.  0  at  X  =  3.5,  5.5  and  7.5. 

(6)  Surface  temperature  vs,  length  at  9  =  30  sec,  and  70  sec. 

(7)  Temperature  distribution  at  x  =  3.5  and  5,5  when  the  surface 
temperature  at  those  points  is  maximum. 

Solution 

Figure  A-5  illustrates  hov;  the  rocket  is  divided  into  regions  for  analysis. 
Further  division  of  the  regions  into  nodes  is  shovra  in  Figure  A-6.  In  general, 
the  entire  length  of  the  rocket  is  divided  into  three  regions  which  are  further 
divided  into  elements  four  inches  long.  Nodes  are  arranged  to  provide  temperatures 
at  a  minimum  of  three  depths  in  the  chamber  wall  and  two  in  the  insulation.  For 
simplification,  only  that  portion  of  the  rocket  shov7n  is  considered  in  the  analysis 
and  all  boundaries  except  those  subjected  to  aerodynamic  heating  is  assumed  to  be 
perfectly  insulated.  In  cases  where  a  more  detailed  analysis  is  required,  the 
entire  structure  plus  propellant  should  also  be  taken  into  account;  hotrever. 
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It  is  not  considered  necessary  in  carrying  out  the  intended  purpose  of  this 
example . 

The  required  thermal  properties  for  Buna-S  and  Spiralloy  rubber  are  taken 

from  Table  III  and  Table  IV,  respectively.  The  contact  coefficient  at  the  chamber 

insulator  interface  is  assumed  to  be  ^  _ Btji _ . 

(sq.ft.)(sec.)(‘’R) 

The  results  obtained  from  the  computer  are  presented  in  the  following  figures 
It  is  obvious  that  the  calculated  temperatures  are  much  too  high  for  proper 
functioning  of  the  rocket  and,  in  a  real  application,  it  would  be  necessary  to 
add  insulation  to  reduce  the  temperature  of  the  critical  components. 

Figure  A-7  Illustrates  the  error  involved  in  using  the  approximate  methods 
of  calculating  heat  transfer  coefficients  for  aerodynamic  heating.  The  range  of 
Mach  numbers  and  altitudes  considered  in  this  example  are  such  that  the  differ¬ 
ences  are  probably  not  critical;  however,  the  margin  of  difference  becomes  larger 
as  Mach  number  and  altitude  are  increased.  It  should  be  noted,  however,  that  the 
resulting  error  in  computing  temperatures  from  the  approximate  film  coefficients 
will  be  smaller,  percentage  wise,  than  the  error  in  the  coefficients. 

Figure  A -8  demonstrates  the  variation  in  film  coefficient  with  time  and 
position.  The  coefficients  are  highest  at  positions  which  are  the  greatest 
distance  from  the  leading  edge.  Another  difference  arises  due  to  the  variation 
in  surface  temperature  along  the  length.  The  reference  temperature  at  which  the 
transport  properties  are  evaluated  for  computing  the  film  coefficient  is  a 
function  of  the  surface  temperature  and  the  boundary  layer  temperature;  therefore, 
the  heat -transfer  coefficient  changes  with  changing  surface  temperature.  The 
boundary  layer  temperature,  shown  in  Figure  A-9,  is  not  a  function  of  axial 
distance. 

Figures  A-10  through  A-12  shew  temperatures  at  selected  times  and  positions 
along  the  length  and  through  the  thickness  of  the  rocket.  Figure  A-13  represents 
the  heat  rate  at  three  positions  along  the  length.  These  figures  contain  only  a 

38 


small  portion  of  the  data  which  were  generated  in  solving  the  example. 

Table  V  shows  the  output  listing  which  is  given  on  the  first  and  last  pages 
of  each  run.  Table  VI  is  an  example  of  the  data  printed  out  at  one  of  the 
thirty-six  times  which  were  called  for  in  this  example. 


PAGE  (FIRST) 


TIMES  AND  MAXIMUM  TEMPERATURES  FOR  TW  NODES 


TABLE  Vi 


TIME 

30.00 

Air. -FT. 

21000.0 

MACH  NO. 

3.000 

PAGE  (  ); 

■ 

Node 

TEMPERATURES-DEG. 

F 

REGION  1 

REGION  2 

REGION  3 

1 

H  2 

.3222E-02 

2. 

0525E-02 

2 

.0203E-02 

1 

+1/3  FT. 

791.8 

SURFACE 

813.8 

-1/3  FT. 

749.4 

,  2 

TW  NODE 

672.1 

808'.  7 

TW  NODE 

626.4 

3 

-1/3  FT. 

8  17*.  7 

803.2 

■Fl/3  FT. 

755.5  ■ 

h 

.075  DEEP 

7  90.7 

798.6 

.075  DEEP 

747.  I 

5 

.075  DEEP 

662-.  8 

TW  NODE 

794.2 

.075  DEEP 

615.5 

6 

.075  DEEP 

817.7 

789.7 

.075  DEEP 

753.7 

- ; - 1 

'  7 

INTERFACE 

789.0 

785.5 

INTERFACE 

744.8 

8 

653.5 

781.5 

60571 

9 

817.5 

777.6 

752.3 

10 

644.3. 

.075  DEEP 

81174 

595 . 2 

11 

817.4 

308.1 

751.9 

12 

INTERFACE 

637.0 

. 

802.5''  - 

'INTERETUrE 

"'587.5  ' 

■  13 

'  INTERFACE 

767.7 

797.7 

INTERFACE 

721.9 

14 

INTERFACE 

602.8 

793.2 

rNTE'RFACE 

5  53. '9 

15 

INNER  'SURF 

764.3 

788.6  , 

INNER  SURF 

718.3 

16 

.  •.'•'I.N.NER  SURF 

597.5 

784.3 

INNE'R  SURF" 

548.8 

17 

.  ■.  ,0.0  NODE 

-459.7  . 

780.2 

0.0  NODE 

-459.7 

.  18 

TBl..,  NODE 

813.0 

776.2 

T'bl  Node 

613.0 

19 

■■  ,NOT  Used 

-459.7 

INTERFACE 

812.1 

NOT  USED 

-459.7 

20 

NOT  USED 

-459.7 

806.8 

NOT  USED 

-■459.7' 

21 

NOT  USED 

-459.7 

801.0 

NOT  USED 

-459.7 

■22 

NOT  USED 

-459.7 

?'96T2 

TrOTTIlED 

=4'5'977 

■2,3 

NOT  USED 

-459.7 

791.7 

NOT  USED 

-459.7 

.24 

N'O'T  USED 

-459.7 

786.9 

'NOr~U'SED'“ 

-459;7 

25 

NOT  USED 

-459.7 

,  , 

782,6 

NOT  USED 

-459,7 

1 

26' 

'NOT  USED 

-'459.7' 

778.5  " 

NOT'  USED  ■ 

-•  rivs'g.T 

/27 

NOT  USED 

-459.7 

774.4 

NOT  USED 

-459.7 

; 

23 

RET  DS'ED 

-'•45977 

■INTERFACE" 

793.2  NOT  USED 

=4^5177 

29 

NOT  USED 

-459.7 

787.5 

NOT  USED 

-459.7 

30 

Tl'dT  USED 

-45577"  ■ 

781.3 

NOT  USED  ■' 

-459.7" 

31 

NOT  USED 

-459.7 

776.2 

NOT  USED 

-459.7 

32 

NOTTrSED 

-459.7 

I7-T-.4 

'NOT  'USED 

-459.7' 

33 

NOT  USED 

-459.7 

766.3 

NOT  USED 

-459.7 

34 

RDTUSE15 

=4!T'9“r7 

/6  1 . 8 

'NOTT.TS'EI! 

-459. 7“ 

35 

NOT  USED 

r459.7 

757.4 

NOT  USED 

-459.7 

36 

NOT  USED 

-459. T 

753.2' 

NOT  USED" 

-459.7 . 

•37 

NOT  USED 

-459.7 

INNER  SURF 

790.2  • 

NOT  USED 

-459.7 

38 

NOT  Used 

-459.7 

784.4 

■ROT  USED 

—459. 7 

'  '  ~  1 

■  39 

NOT  USED 

-459.7 

778.1 

NOT  USED 

-459.7 

’  -40 

NOT  USED 

-459.7 

77170 

"NOTUS'E'D 

-459 . 7 

41 

NOT  USED 

-459.7 

768.1 

NOT  USED 

-459.7 

42 

NOT  USED 

'  -459. 7~ 

■■  763.0  " 

NOT  USED  " 

■ -45'9.7 

43 

NOT  USED 

-459.7 

758.4 

NOT  USED 

-459.7 

44 

NOT  USED 

-459.7 

'  754.0'. 

NOT  USED  ' 

■  -'459.7  ■ 

4-5 

NOT  USED 

-459.7 

749.8 

NOT  USED 

-459.7 

46 

NOT  USED 

-459.7 

0.0  NODE 

-4'59.'7 

UO’T'USE'D 

-45977 

47 

' NOT  USED 

-459.7 

TBL  NODE 

813.0 

NOT  USED 

-459.7 

TBL,.TW,Q 

RATE  18  2 

3.27 

a't  5 

0.38 

18  2 

3.77 

EXAMPLE  III 


A3  POLARIS  Second-Stage  Motor  Surface  Temperature 
Due  to  Aerodynamic  Heating _ _ 

The  outside  motor  case  temperature  of  the  A3  POLARIS  second-stage  motor  was 
calculated  as  a  function  of  time  at  missile  station  220  using  A3X-2  flight 
trajectory  data.  These  temperatures  were  compared  with  measured  temperatures 
obtained  at  the  same  station  during  the  A3X-2  flight.  The  measured  outside 
motor  case  temperature  data  are  compared  with  calculated  values  in  Figure  A-14. 
Due  to  the  variation  in  reported  values  of  thermal  conductivity  for  Spiralloy, 
the  calculations  were  carried  out  using  two  different  values;  results  for  both 
are  presented  for  comparison  in  Figure  A-14.  It  will  be  noted  that  reasonable 
agreement  between  the  actual  flight  data  and  the  calculated  values  vaa  obtained. 

The  sudden  temperature  rise  in  the  actual  case  temperature  at  a  time  of 
one  second  was  caused  by  first-stage  ignition  and  should  not  have  had  an  effect 
on  the  peak  temperature  measured  at  56  seconds. 


TIME  -  HOURS'  .2893 

FIGURE  A1  ^ 

Surface  Temperature  History 


43 


Note !  Points  correspond  to 

‘"Analytical''  Solution  taken 
from  Schack  Charts  using 
Average  Properties. 


FIGURE  A4 

Schemaflc  Representafion  of  fhe  Propulsion  Unif  for  a  Hypofheticof- Afr^To-Air  Missile 
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•  FIGURE  A5 

Division  of  Rocket  Into  Regions  for  Aerodynamic  Heating  Analysis 
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Region  "G"  Subdivision  is  a  mirror  image  of  Region  "A”.. 
Nodes  1,4, 7-, 13  and  15  in  Region  "C  connect  with  nodes  9,18, 
27,36,  and  45,  respectively,  in  Region  "B". 


FIGURE  A6 

Division  of  Regions  info  Nodes  for  Aerodynamic  Keafing  Analysis 
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FIGURE  A8 

Coefficient’  at  Three  Positions  Along  the  Axis  of  a  Hypothetical  Missile 


FIGURE  A9 

Boundary  Layer  Temperature  fora  Hypothetical  Missile 


FIGURE  AlO 

Propel  lant“~  IjssssI^flsDr  Inteifece  Temperatures 


AXIAL  DISTANCE  FROM  THE  LEADING  EDGE 


